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Abstract. Dawn space mission will provide the first, detailed data of two of the 
major bodies in the Main Belt, Vesta and Ceres. In the framework of our studies 
on the origin of Solar System, we modelled the accretion of Jupiter and, through 
an N-Body code developed on purpose, we evaluated the flux of impactors on Vesta 
and Ceres keeping track of their formation zones. We also studied the effects of 
the possible inward migration of Jupiter on the rate and the characteristics of the 
impacts. We here describe the different scenarios and their implications for the 
evolution of Solar System. 



1. Introduction 

Dawn mission will provide detailed images of Vesta and Ceres surfaces and sup- 
ply crucial information to constrain their mineralogical and elemental composition 
through VIR, its imaging spectrometer. Thanks to these data, we will be able to 
study in depth the crater record on the surface of both Vesta and Ceres. We discuss 
here our study of the cratering process of Vesta and Ceres at the time of Jupiter 
formation. In our model we consider Jupiter's gas accretion and displacement due 
to angular momentum exchange with the surrounding n ebula. The gas accretion 
model used has been described in Coradini et al. ] (12004 ) . Jupiter's migration has 



been included to estimate the effect of increasing displacements. 

2. Dynamical and physical model 

To explore the early collisional history of Vesta and Ceres we simulated the dy- 
namical evolution of a section of the young Solar System at the time of Jupiter's 
core formation and the subsequent accretion of the gaseous envelope. Our template 
of the forming Solar System was composed of the Sun, the accreting Jupiter and a 
swarm of massless particles representing the planetesimals. The massless particles 
were initially distributed into a limited spatial region, which has been chosen to op- 
timise the computational requirements after a set of numerical experiments aiming 
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to determine the region of Solar System influenced by the forming Jupiter on the 
considered timespan. Being interested on the effects of Jupiter's mass increase on the 
dynamical stability and the evolution of the Solar System, we didn't model the giant 
planet formation process directly through hydrodynamical computations. In our sim- 
ulations we reproduced the evolution of Jupiter through an analytical approach: the 
parameters on which the model was based, however, were derived fr om the results of 



hydro dynamical simulations performed with the code described in ICoradini et al. 



(12004 ) . During the dynamical evolution of our template of the Solar System we 
evaluated through a statistic approach the probability of planetesimals impacting 
against Vesta and Ceres. In the following subsections we will describe in detail the 
initial conditions and the physical parameters and constrains of the model. 

2.1. Jupiter's formation and migration. In our analytical approach, we consid- 
ered Jupiter's formation as composed by 3 different stages: 

• a first core accretion phase; 

• a second exponential gas accretion phase; 

• a final phase of asymptotic gas accretion. 

At the beginning of the simulations, Jupiter is an embryo with mass Mo = 0.1 M® 
growing on a time r c = 3 x 10 5 years to the critical mass M c = 15 M®. The mass 
growth is governed by the equation 

(1) M p = M + (-^j (M c - M ) x (l - e */^) 

The timescale for core accretion, smaller than what the theoretical models for giant 
planet formation predict, has been chosen to limit the computational load of the 
simulations. During the first stage of accretion, Jupiter's core mass is negligible in 
terms of its perturbing effects on the planetesimals in the disk and the only plan- 
etesimals affected by Jupiter are those which undergo a close encounter with the 
forming planet. As a consequence, in this phase the number of impacts against the 
considered asteroids is governed by stochastic collisions and is a linear function of 
time: we can approximate the evolution on a longer timescale by extrapolating from 
our results. 

After Jupiter's core reached the critical mass value of 15 M®, the nebular gas sur- 
rounding it becomes gravitationally unstable and starts to be rapidly accreted by the 
planet to form its massive envelope. During this phase, Jupiter strongly perturbs 
near-by planetesimals, clearing a gap of increasing width in the planetesimals' disk. 
To model this second phase of Jupiter's form ation we used a t i mesca le derived from 
the hydrodynamical simulations described in Coradini et al. ( 20041 ): the value we 



employed in the simulation is r g = 5 x 10 years and is linked to the turbolence pa- 
rameter a = 0.01 emploied in modeling the disc. During this phase, Jupiter's mass 
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growth is governed by the equation 



(2) 



M p = M c + (Mj - M c ) x (l - e (i - Tc)/T ») 



where Mj = 1.8986 x 10 30 ^ = 317.83 M ffi is Jupiter's final mass. 
In the third and final stage, Jupiter's accretion slows down while the giant planet 
reaches its final mass value and its gravitational perturbations start to secularly affect 
more and more distant planetesimals. While Jupiter's mass growth is still governed 
by eqj2] and there is no discontinuity from a numerical point of view between this 
phase and the previous one, we considered them as two separate stages of Jupiter's 
formation due to their different effects on the early Solar System, as we will describe 
broadly in section [3J We followed the evolution of the system across the second and 
third stages for r a = 7 x 10 5 years. Theoretical models indicate that forming giant 
planets go subject to external torques caused by the gravitational interaction between 
the planet and the protoplanetary disk: as a consequence of these torques, Jupiter is 
expected to radially migrate inward. To evaluate the effects of the radial displacement 
of the planet on the dynamical stability of the planetesimals' disk, in our simulations 
we forced Jupiter to migrate. Like the growth of the gaseous envelope, also the radial 
displacement starts at the beginning second phase of Jupiter's formation and reaches 
an asymptotic state in the third phase. The equation governing the radial migration 
is similar to eq. [2J 



where r$ is Jupiter's position at the beginning of the simulation (depending on the 
desired displacement), rj is the final (present) position and r r = r g = 5 x 10 3 years. 
In order to better evaluate the effects of migration, in our set of simulations we 
considered four different displacements: AU (no migration), 0.25 AU, 0.5 AU and 
1 AU. The initial semimajor axis of Jupiter has been consequently chosen so that its 
value at the end of the simulations would correspond to the one the planet possesses 
at present time. 

2.2. Dynamical characterisation of the planetesimals. The initial distribution 
of the planetesimals in our simulations populated a region of phase space defined as 
follows: 



(3) 




(4) 



2 AU < ai < 10 AU 

< a < 3 x 10~ 2 
rad < %i < 10~ 2 rad 
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The values of eccentricity and inclination associated to each massless particle in our 
simulation have chosen randomly as 

(5) ti = e X 
and 

(6) U = io(l-X) 

where eo = 10~ 2 , Zq = 10~ 2 rad and X is a random number extracted from a uniform 
distribution with < X < 1. The radial distribution of the planetesimals has 
been obtained by assuming the relations hip between the mass and the heliocentric 



distance of the planetesimals described in ICoradini et al. I (119811 ) with superimposed 



a Maxwell-Boltzmann distribution. The number surface density in the planetesimals' 
disk can therefore be expressed as a function of mass and radial distance as 

(7) n(m, r) = Q(r)mV (m/ ™ p(r))2 

where Q(r) re presents the rad i al dep endence of the number surface density n(m,r) 



3P_r 

and, following ICoradini et al. I (1l98ll ). m p (r) is the average mass of the planetesimals 



as a function of heliocentric distance given by the relationship 



m p = m 



IAU 



where m p and mo are expressed in g, r is expressed in AU and (3 = 1.68. The value 
m is the average mass of a planetesimal at IAU. The functional form of Q(r) can 
be obtained by coupling eq. [7J with the relationship governing the surface density 
profile in the Solar Nebula 

(9) Gpir) — cta{r) = aao 



IAU 



where we assumed the standard density-heliocentric distance relationship and put 
n s = 1.5. The value a p is the mass density of the planetesimals, a is the mass density 
of the gas in the Solar Nebula, <To is the mass density of the gas at 1 AU and a is the 
dust to gas ratio: here we made the classical assumption a = 10~ 2 . In the following, 
all the distances r, when normalised to 1 AU, will be indicated with the capital letter 
R while the symbol AU will indicate the value IAU = 1.49597870691 x 10 13 cm. 
By integrating eq. [7] over m we obtain 



(10) n*(r) = / n(m,r)dm = 

Jo 

= J™ Q(r)m 2 e- {m/m ^ r))2 dm = ^Q(r)m p (r) 
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while for the mass surface density of planetesimals we have 
(11) Op(r)= / n(m,r)mdm= I Q(r)mV (m/mp(r)) dm = -Q(r)rr? v (r) 



o Jo 2 

By equating eq. [11] with eq. [9] we get 

(12) Q(r) = 2^iT"« 

which, substituting [12] into [10] and applying eq. gives 



(13) n*(r) = xJL^l R -in s+ 0) 

2 mo 

Now we can obtain the cumulative distribution of n*(r) by integrating 

(14) N(x) = f X 2iTrn*(r)dr 



3/2« ff (ATT ,2 ( 1 \ /V X \ 2-n s -/3 , r x 2-»,-/J\ 

where r m j„ = 2 At/ is the inner radius of the planetesimals' disk. We can obtain the 
total number of planetesimals in the disk by integrating eq. [H] between r m j n = 2 AU 
and r max = 10 At/ : 



(15) N tot = / 2irrn*(r)dr = 

= vr 3 / 2 ^ (AU) 2 ( 1 ^ iV^V"^ _ f^V~ ns ^ 
m 1 ; \2-n s -f3) \\ AU J \ AU J J 

The function X = uniformly varies in the range < X < 1. Once we invert 
previous equation and express the initial position of the planetesimals in the disk as 
x = f(X) where 

n s +/3~2 

(16, *=fm=[x{{^y~'-{%y "x^ 1 " "" 

we can use X as the random variable in a Montecarlo extraction to spatially populate 
the disk. Through eq. Owe can link the average mass m p {r) of the planetesimals to 
their initial position r = x (which we assume coinciding with the formation region). 
To obtain each planetesimal's mass value, we apply again a Montecarlo method 
through a uniform random variable Y. To obtain Y we need to compute and nor- 
malise the cumulative distribution of eq. [7] considering m as the variable of interest 
and T clS db constant. The cumulative distribution is 

(17) r n(m,r)dm=^^-R- n ° ( V y 2 e- y2 dy 

Jo mJr) Jo 
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where y = m/m p {r) and y* = m*/m p (r). 

The integral in the right side of eq. [T7] cannot be solved analytically: substituting 
z = y 2 with dz = lydy we obtain 



(18) \ V \\-y'dy= \ l'~ z l ' 2 e~ z dz 

Jo I Jo 

The integral on the right side of previous equation is a lower incomplete Gamma 
function with real parameter a = 3/2 or 7 ^3/2, \/z*). Normalising eq. [T71 over eq. 
[I3]we get the uniform random variable Y varying in the range [0, 1] 

(19) y ^ 2 7 (3/2,y*) ^ p 

'7T 



where P(3/2,y*) is the lower incomplete Gamma ratio. The inverse of the lower 
incomplete Gamma ratio can be computed numerically and, by substituting y* back 
with m* /m p (r) we obtain 

(20) m(r) = rfipinv (P (3/2, Y))) 

We can therefore assign to each planetesimal its own mass value through a second 
Montecarlo extraction. 

The number of massless particles we used in our simulations has been n mp = 8 x 10 4 to 
optimise the computational load. The use of massless particles assures the linearity 
of the processes over the number of considered bodies: we can therefore extrapolate 
the number of impacts expected in a disk of planetesimals described by the density 
profiles previously presented by multiplying the number of impacts recorded in our 
simulations by a factor 7 where 

(21) 7 = N tot /n mp 

where N to t is computed through eq. [151 The 7 factor is independent from heliocentric 
distance r since the n mp planetesimals have been distributed following the number 
density profile n*(r) of the disk. 

The dynamical evolution of Jupiter and the swarm of massless particles is computed 
through a fourth order Runge-Kutta integrator with a self adjusting time-step. To 
choose the time-step, our numerical code evaluates the smallest timescale in the 
simulations at a given time between: 

• the orbital periods of the massless particles 

• the orbital periods of Jupiter and the two asteroids 

• the free-fall timescale of each pair massless particle- Jupiter considered as an 
isolated system 

The time-step employed is computed by the code as 

(22) ti s T m i n / ft s 
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where f ts = 100 in our simulations. 

2.3. Chemical and physical characterisation of planetesimals. To evaluate 
the possible amount of volatiles delivered to the asteroids from the planetesimals, 
we recorded the initial positions of the massless particles in our simulations. We 
assumed that the initial positions of the particles coincided with their formation 
regions. We considered two compositional classes: 

• particles formed inside the Snow Line (ISL in the following) were considered 
rocky bodies 

• particles formed beyond the Snow Line (BSL in the following) were consid- 
ered icy bodies 

ISL and BSL particles were characterised by mean density values respectively pisl = 
3.0 g/cm 3 and Pbsl = 1.0 g/ cm 3 . T he location of the Snow Line in our simulations 
has been fixed at tsl = 4.0 AU (see Encrenaz ( 20081 ) and references therein). 



2.4. Collisional history. To reproduce the collisional histories of the two asteroids, 
instead to follow the explicit dynamics of Vesta and Ceres and proceed in a deter- 
minist way, we opted for a statistical approach. Each asteroid was spread into a 
torus with mean radius R? and section ot defined as 

(23) R T = a A 
and 

(24) a T = j x (D A f G f 

where a A and D A were respectively the semimajor axis and the physical diameter of 
the considered asteroid while fa is the gravitational focusing factor 

(25) f a = l+ ' 

\ ^enc 

with v esc being the escape velocity from the asteroid and v enc being the relative 
velocity between the asteroid and the impacting planetesimal. The gravitational 
focussing factor fa was introduced to account for the perturbations produced by the 
asteroid on the orbits of the planetesimals, which were not included in the explicit 
dynamical model. 

When a planetesimal crosses one of the asteroids' torii, the impact probability is 
the probability that both planetesimal and asteroid will occupy the same region at 
the same time. This probability can be evaluated as the ratio between the effective 
collisional time and the orbital period of the asteroid. The effective collisional time is 
the amount of time available for collisions and is evaluated as the minimum between 
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the time spent by the asteroid and the planetesimal into the crossed region of the 
torus. This is equivalent to writing 

(26) P„„ = mm ^ TA) 

where Ta is the orbital period of the asteroid, ta is the time spent by the asteroid 
into the crossed region, Tp is the time spent by the planetesimal into the torus and 
min{jp, ta) is the effective collisional time. 

To estimate r P and ta we need to identify the intersections of the planetesimals' 
orbits with the asteroids' torii. Our approach is based on a linearisation of the 
problem, i.e. we approximated with a line the orbital path of the planetesimal while 
passing near a torus and solved the ray-torus intersection problem. A torus centred 
on the origin of the axes and lying on the xy plane is described by the equation 



2 2 
Z — S 



(27) f(x, y, z) = (x 2 + y 2 + z 2 - (R 2 + s 2 )) + 4R 2 ( 

where R = Rt and s = 0.5 (DaJg)- The surface of the torus is described by the 
homogeneous equation f(x, y, z) = 0: by substituting the variables (x, y, z) with the 
vectorial equation of the ray associated to the linearised orbit 

(28) r = rht + q 

where m and q are respectively the angular coefficient and the origin of the ray, we 
obtain the fourth order equation identifying the intersections between the ray and 
the torus. By writing 

(29) a = rh ■ rh 

(30) b = m-q 

(31) c = q-q 
we can express the fourth order equation as 

(32) rf 4 t 4 + d 3 t 3 + d 2 t 2 + dxt 1 + d = 
where 



(33) 



d\ 


= a 2 




d 3 


= 4ab 




d 2 


= 4& 2_ 2as 2_ 


2aR 2 + lac + 4R 2 m 2 z 


di 


= 46c -46s 2 - 


AbR 2 + 8m z q z R 2 


do 


= c 2 + R 4 + s 4 


- 2R 2 s 2 - 2cs 2 - 2cR 2 + 4R 2 q 2 
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This fourth order equation can be solved analytically through Ferrari's method. The 
quartic equation will have up to four real solutions, depending on the number of 
intersections: by selecting the appropriate pair of real solutions, if existing, we can 
derive the coordinates of the intersection points x\ and x%. From the intersection 
points we can derive the length of the planetesimal's path d P through the torus and 
the crossing time Tp which, in the linear approximation, is 

(34) tp = ^ 

\Vp\ 

where \vp\ is the modulus of the velocity of the planetesimal. 

From the knowledge of the intersection points we can also derive the angular width 
of the crossed section of the torus AO a and the time ta the asteroid spends into the 
crossed region: 

T A A6 a 



(35) t a = AO 



A 



2n uja 



where uja = Ua = yt is the angular velocity (coinciding with the orbital mean motion 
for circular orbits) of the asteroid. 

2.5. Impact cratering. To evaluate the effects of the impacts on the surface crater- 
ing of the two asteroids, for each recorded collision we computed the averag e crate r 
diameter using the empirical scaling law (see p. 165 of De Pater fc Lissauer ( 200ll )) 

/p X 0.22 

(36) D = 1.8p°- n P y 3 (2R*f 13 ( — J (stnO) 1/3 

where the letters a and i indicate the asteroid and the impactor respectively, p is 
the density, R the physical radius, g the gravitational acceleration, is the impact 
kinetic energy and is the impact angle respect to the local horizontal. In our 
estimation of the surface cratering w e used a fixes value = 45° as the average 



impact angle (see iPierazzo fc Melosh I (120001 ) and references therein) . 



2.6. Parallelising the model. Before concluding the section devoted to the dy- 
namical and physical model, we would like to emphasize one advantage of our mod- 
elization, i.e. the possibility to parallelise the evolution of the template of the Solar 
System in a straightforward way. Since the disc of planetesimal is composed only of 
massless particles, which by definition do not interact between themselves, and since 
Jupiter's evolution is driven by the migration and accretions rates described by eqs. 
(U E] and [3l we are able to split each simulation into a set of sub-problems which 
can be treated in parallel. In our implementation, we used a GRID-like approach 
by dividing the disc of planetesimals into a number of concentric rings containing a 
fixed amount of test particles. Each ring has been evolved independently under the 



10 



TURRINI D., MAGNI G., AND CORADINI A. 



influence of the forming Jupiter and the collisions with the two asteorids recorded. 
At the end of the simulations, by merging the results we obtained a representation 
of the evolution of the system as a whole. Obviously, a different normalisation factor 
has to be used for each ring of massless particles to normalise it to the real popu- 
lation expected in the orbital region considered. Through this approach, we have 
been able to run the equivalent of a 3- month long simulation with 8 x 10 4 massless 
particles by running a set of 8 sub-simulations with 10 4 massless particles each, every 
sub-simulation taking about 12 days to conclude. 

3. Results 

In our simulations we identified 4 distinct dynamical families contributing to the 
cratering of Vesta and Ceres in the temporal framework considered. Two of these 
families (a resonant and a non-resonant one) belonged to the ISL population and 
two (again, a resonant and a non-resonant one) to the BSL population. 

3.1. Dynamical features of the impactors. The ISL impactors can be divided 
into two main dynamical groups (see figs. [I] and [2]). The first one is composed by 
bodies orbitally near the two asteroids and represents the population of bodies con- 
tributing to the accretion process of the asteroids, hereafter indicated as primordial 
impactors. For this population collisions start shortly after the beginning of the 
simulations and would deplete the population if not for the dynamical perturbations 
of Jupiter, which inject new massless particles on Vesta-crossing or Ceres- crossing 
orbits. The second ISL group is composed of families of resonant impactors, with 
the number of active resonances and their relative abundances of impactors strongly 
depending on the extent of Jupiter's displacement. As a rule of the thumb, the res- 
onance affecting each asteroid for all cases considered of Jupiter's migration is the 
one located nearest to the respective orbital regions. The flux of impactors due to 
the resonances increases with increasing displacements of Jupiter (see figs. [1] and 
[2]), as is to be expected since the resonances co-move inward with the giant planet 
and sweep larger orbital regions, affecting more massless particles. The resonances 
affecting Vesta and Ceres in the temporal framework of our simulations are located 
at about 2.5 AU and 3.5 AU (see figs. [T]and[2]): they appear in the simulations after 
the beginning of the rapid gas accretion phase, when Jupiter's mass is high enough to 
strongly excite the resonances. Ceres' resonant impactors (fig. [2]) have eccentricities 
distributed between 0.1 — 0.2 if originating from the 2.5 AU region and 0.1 — 0.6 
if originating from the 3.5 AU region. The eccentricity values of Vesta's resonant 
impactors (fig. [1]) are distributed in the range 0.05 — 0.2 for objects coming from the 
2.5 AU region and 0.2 — 0.5 for those coming from the 3.5 AU region. 
Also between the BSL impactors we detected two dynamical groups (see figs. [T] 
and [2]): a first, sparse one composed of non- resonant bodies scattered by Jupiter's 
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Figure 1. Semilogarithmic plots of the orbital elements of Vesta's 
impactors in the a — e plane. Red and green symbols represent respec- 
tively the final (at impact) and initial orbits of ISL impactors, while 
the blue and magenta symbols represent respectively the final and 
initial orbits of BSL impactors. As is clearly visible from the plots, 
ISL impactors are divided in two dynamical families, a first, primor- 
dial one and a second, resonant one (see section [37T1 for details). Also 
the BSL impactors can be divided between non-resonant and resonant 
bodies, but the latters disappear for increasing values of Jupiter's ra- 
dial displacement (0.5 and 1 AU). Note that these plots show only the 
dynamical classes of impactors recorded in the simulations since they 
are not normalised to the real disc population. 



gravitational perturbations and a second, narrower one (located near the Snow Line) 
related to a mean motion resonance with the giant planet. Non-resonant planetesi- 
mals are dynamically destabilised when Jupiter's mass reaches a significant fraction 
of its final value and part of them is injected into the inner Solar System on Ceres- 
crossing and Vesta-crossing orbits. The existence and the abundance of the BSL 
resonant population are strictly linked to the extent of Jupiter's radial displacement, 
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Figure 2. Semilogarithmic plot of the orbital elements of Ceres' irn- 
pactors in the a— e plane. Red and green symbols represent respectively 
the final (at impact) and initial orbits of ISL impactors, while the blue 
and magenta symbols represent respectively the final and initial orbits 
of BSL impactors. As in fig. [TJ both ISL and BSL impactors can 
be divided into a non-resonant and a resonant population. In Ceres' 
case, however, the resonant BSL population is present in three of our 
dynamical scenarios (i.e. it doesn't appear in the case of Jupiter's 
0.5 AU radial displacement). Note that these plots show only the dy- 
namical classes of impactors recorded in the simulations since they are 
not normalised to the real disc population. 



both for Vesta and Ceres. Increasing values of the radial displacement cause a de- 
crease in the number of resonant BSL impactors. However, in Ceres' case, a new 
resonant population originating near 4.6 AU appears in the 1 AU migration scenario. 
The final orbits of all BSL impactors are characterised by high eccentricity values, 
in the range 0.3 — 0.7 for Vesta (fig. [TJ and 0.2 — 0.6 for Ceres (fig. [2J). 
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Figure 3. Frequency versus formation region histograms of Vesta's 
impactors normalised over the real planetesimal population: red bars 
are those related to ISL impactors, blue ones to BSL impactors. The 
frequency of ISL impactors in all plots has been divided by a factor 100 
to enhance the readability and facilitate the comparison of the plots. 
Once normalised, ISL planetesimals completely dominate the flux of 
impactors on Vesta. The peak in BSL impactors in the top left plot is 
due to the excitation of an orbital resonance with Jupiter. 



3.2. Characterization of the flux of impactors. The dynamical classes of im- 
pactors we showed in figs. [T] and [2] and discussed in previous section represent the 
events recorded during our simulations and are not representative of the real flux 
of impactors. In figs. [3] and H] we show the normalised frequency versus formation 
region distributions of both ISL and BSL impactors in all the migration scenario we 
considered. As is clearly visible, ISL impactors are 2 — 3 orders of magnitude more 
abundant than BSL impactors and completely dominate the craterisation histories 
of the two asteroids. The frequency of BSL impactors is greatly enhanced in those 
migration scenarios where Jupiter excites resonances in the orbital region outside the 
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Figure 4. Frequency versus formation region histograms of Ceres's 
impactors normalised over the real planetesimal population: red bars 
are those related to ISL impactors, blue ones to BSL impactors. The 
frequency of ISL impactors in all plots has been divided by a factor 
100 to enhance the readability and facilitate the comparison of the 
plots. The flux of BSL impactors, while by far inferior than that 
of ISL bodies, is significantly higher than in the case of Vesta. The 
peaks in BSL impactors in the top left and bottom right plots are due 
to the excitation of different orbital resonances with Jupiter, while 
the importance of the analogue resonance in the 0.25 AU migration 
scenario is greatly reduced after the normalisation. 



Snow Line, as we mentioned in previous section. 



3.2.1. ISL impactors. The flux of ISL impactors on the two asteroids is characterised 
by two phases (see figs. [5] and [6]). The first phase covers the time it takes the planet 
to form a critical mass core and accrete a significant mass of gas: during this phase 
the flux of impactors on the asteroids is dominated by primordial impactors. The 
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Figure 5. Normalised cumulative distribution of the flux of impactors 
over time for Vesta: the red curve is related to ISL impactors, the blue 
one to BSL impactors. In all plots the fluxes have been normalised to 
their values at the end of the simulations. The flux of ISL impactors 
is divided into two main phases: a first one due to the contribution of 
primordial impactors and a second one, starting after Jupiter accretes 
a significant mass, due to both primordial and resonant impactors. 
In both phases, the growth of the cumulative flux slows down with 
time and eventually stops. The flux of BSL impactors begin when 
Jupiter's gaseous envelope reaches a significan mass and the cumulative 
distribution grows until the end of the simulations. 



steepness of the cumulative flux during this phase varies depending on Jupiter's 
location: the farther Jupiter formed from its present orbit, the lesser the steepness of 
the curve (see figs. [5] and [6]). This can be intuitively explained by Jupiter being less 
effective at injecting new primordial impactors on Vesta-crossing or Ceres- crossing 
orbits when it forms farther from their location. In the last third of the first phase, 
between 2 x 10 5 and 3 x 10 5 years, the flux of ISL primordial impactors stops both 
in Vesta's and Ceres' cases (see figs. and El). The lenght of the stop in the flux of 
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Figure 6. Normalised cumulative distribution of the flux of impactors 
over time for Ceres: the red curve is related to ISL impactors, the blue 
one to BSL impactors. In all plots the fluxes have been normalised to 
their values at the end of the simulations. Most of the considerations 
exposed for Vesta's case (see fig. [5]) are valid also for Ceres. The tran- 
sition between the two phases for ISL impactors mentioned in fig. [5] is 
clearly dependent from Jupiter's formation region. For Jupiter forming 
near its final semimajor axis (0 and 0.25 ^£7 cases), the gravitational 
perturbations of the growing core can affect the ISL planetesimals and 
cause a visible enhancement in the flux of impactors. If Jupiter formed 
farther (0.5 and AU cases) the transition between the two phases is 
more clean. 



impactors depends on the initial position of Jupiter (see figs. [5] and [6]). If Jupiter 
formed nearer to the asteroids, as in the and 0.25 AU migration scenarios, the core 
would start to inject new bodies on asteroid-crossing orbits once it approaches the 
critical mass. The influence of Jupiter's core then results in a smoother transition 
between the first and the second phase. If Jupiter formed farther away, as in the 0.5 
and 1 AU migration scenarios, the core is unable to accomplish the same task and 
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Figure 7. Cumulative distribution of the mass flux over time for 
Vesta: the red curve is related to ISL impactors, the blue one to BSL 
impactors. The mass flux of ISL impactors in all plots has been divided 
by a factor 100 to enhance the readability and facilitate the comparison 
of the plots. The overall features of the cumulative distribution of the 
mass flux are similar to those of the cumulative distribution of the flux 
of impactors (see fig. [5] for details). 



the stop lasts till the beginning of the second phase, resulting in an more abrupt 
transition. The second phase starts when Jupiter is massive enough to significantly 
perturb the orbital resonances (see figs. OandE]). The flux of ISL impactors greatly 
increases on timescales of the order of that of Jupiter's mass growth, almost doubling 
in about 2 — 3 x 10 years. The flux then slows down as Jupiter clears the reservoir of 
planenesimals located near the mean motion resonances and bodies farther away are 
perturbed on longer timescales. During this second phase, the flux of ISL impactor 
is still sensitive to Jupiter's migration. This dependence can be seen from figs. [7] and 
[HJ the figures show the cumulative mass fluxes on Vesta and Ceres respectively and 
the final values of the mass flux increase linearly with Jupiter's radial displacement. 
Finally, even if similar when normalised (see the captions of figs. and for details), 
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Figure 8. Cumulative distribution of the mass flux over time for 
Ceres: the red curve is related to ISL impactors, the blue one to BSL 
impactors. The number of ISL impactors in all plots has been divided 
by a factor 100 to enhance the readability and facilitate the comparison 
of the plots. The overall features of the cumulative distribution of the 
mass flux are similar to those of the cumulative distribution of the flux 
of impactors (see fig. [6] for details). 



the flux of impactors (both primordial and resonant) on Vesta is about half than the 
corresponding one on Ceres (see figs. [7] and [5]) • This is due to the inner position of 
the former asteroid, which imposes stricter constrains on the dynamical evolution of 
the planetesimals to be injected on Vesta-crossing orbits. 

3.2.2. BSL impactors. The flux of BSL impactors greatly varies depending on the 
considered asteroid and migration scenario. For Vesta, the only relevant contribution 
of BSL impactors to the cratering history of the asteroid appears in the AU, no- 
migration scenario (see fig. [3]). In this scenario, Vesta is hitted by about 2 x 10 2 
impactors (see fig. [3]) once the flux is normalised on the real population of the 
disk of planetesimals. In all the other scenarios, the contribution of BSL impactors 
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is of the order of 10 events (see fig. [3]). This implies that the expected flux of 
volatile-rich bodies on Vesta at the time of Jupiter's formation is about 3 orders of 
magnitude lower than that of rocky bodies, with the only exception being the case 
where Jupiter formed at its present location. Due to the mass-heliocentric distance 
relationship assumed in eq. [HI the mass of the planetesimals in the disc increases with 
the heliocentric distance. Once this effect is taken into account, in the no-migration 
scenario the mass flux of BSL planetesimals is about a factor 40 lower than the one 
of ISL planetesimals (see fig. [7]). In all other three migration scenario, the ratio 
between the mass fluxes of ISL and BSL planetesimals doesn't change significantly 
(see fig. ED. 

Ceres' case is analogous to Vesta's one, but the BSL flux of impactors is enhanced in 
the AU and in the 1 AU cases due to the excitation of mean motion resonances. The 
enhancement results in overall fluxes 1 order of magnitude and mass fluxes 2 orders 
of magnitude greater than in the 0.25 AU and the 0.5 AU cases (see figs. H]and[B]). 
The enhanced fluxes and mass fluxes of BSL impactors are about a factor 3 higher 
than in Vesta's enhanced case (OAU). In the non-enhanced cases (i.e. 0.25 AU and 
0.5 AU) the mass fluxes of BSL impactors on Ceres are about an order of magnitude 
higher than the ones in Vesta's non-enhanced cases (see figs. [7] and [8]). 
The overall delivery of volatile-rich material in the orbital region of Ceres remains 
more than one order of magnitude lower than that of rocky material. We must note, 
however, that the injection of BSL planetesimals on asteroid-crossing orbits, being 
a secular process, didn't stop before the end of the simulations (see the trend of the 
BSL curves in figs. [5] and ED- The secular perturbations of Jupiter (and, in the real 
Solar System, of the other giant planets) could inject more bodies on asteroid-crossing 
orbits on timescales longer than the one covered by our simulations. 

3.3. Characterization of the impacts. In our simulations we recorded the char- 
acteristics of each occurred impact detailing them by their dynamical features and 
their effects on the target (i.e. crater size and energy delivered to the asteroid). In 
the following we will describe our results. 

3.3.1. Dynamical features. While considering the dynamical features, we observe 
three classes of impacts which are common both to Vesta's and Ceres' impactors. 
The first class is composed by the ISL primordial impactors: is characterised by 
low impact velocities, always inferior to 1 km/s, and a continuous mass spectrum 
over the range we considered (see figs. [9] and [10]). The second class is composed by 
the ISL resonant impactors: this class too has a continuous mass spectrum, but the 
characteristic velocities are systematically higher than those of the first class, rang- 
ing between 1 — 10 km/s (see figs. 191 and [TO"]) . The third class of impacts is due to 
the BSL impactors: there is no evident distinction between the impacts due to BSL 
resonant and non-resonant planetesimals, with the only exception of the Ceres' 1 AU 
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Figure 9. Dynamical characterisation of the impacts on Vesta: red 
symbols indicated ISL impactors and blue ones BSL impactors. The 
data reproduced in the plots are the events recorded in our simulations: 
they have not been normalised to the real population of the planetes- 
imal disc and are used only to show the observed dynamical classes. 
As can be seen from the plots, ISL resonant and primordial impactors 
are characterised by different dynamical features, while no difference 
is evident between resonant and non-resonant BSL impactors in those 
scenario with Jupiter migrating. In the no-migration scenario BSL im- 
pactors seem to divide in two clusters, a first between 3 — 4 km/ s and 
a second between 8 — 10 km/ s. The relative abundance of bodies in 
the two clusters, however, argues against attributing this effect to the 
non-resonant versus resonant subdivision (see fig. HJ 



migration scenario. The BSL impactors globally tend to cluster towards higher mass 
values as a consequence of the mass-heliocentric distance relationship assumed in the 
model (eq. |SJ). The impact speeds of BSL impactors mostly overlap the same range 
of values of the resonant ISL impactors but, depending on the migration scenario, 
a fraction can hit the asteroids with velocities higher than 10 km/s (see figs. Oand 
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Figure 10. Dynamical characterisation of the impacts on Ceres: red 
symbols indicated ISL impactors and blue ones BSL impactors. The 
data reproduced in the plots are the events recorded in our simulations: 
they have not been normalised to the real population of the planetes- 
imal disc and are used only to show the observed dynamical classes. 
As can be seen from the plots, ISL resonant and primordial impactors 
are characterised by different dynamical features, while no difference is 
evident between BSL resonant and non-resonant impactors in the first 
three migration scenarios. In the 1 AU migration scenario, resonant 
BSL impactors cluster in the range 2 — 3 km/ s while non-resonant BSL 
impactors are characterised by higher velocities. 



ITOl . In one case we recorded an impact velocity value of about 20 km/s (see fig. [TU1 
bottom right panel). 



3.3.2. Impact craters. While considering the craterisation due to the impacts, we 
observe that ISL and BSL planetesimals form two different classes. As a general 
rule, ISL-related craters completely dominate the crater distributions on the two 
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Figure 11. Normalised frequency versus crater size histograms of 
Vesta's impactors: red bars are those related to ISL impactors, blue 
ones to BSL impactors. The frequency of ISL impactors in all plots 
has been divided by a factor 100 to enhance the readability of the 
plots. As can be seen: the high-end tail of the ISL crater size distribu- 
tion moves towards bigger sizes for increasing values of Jupiter's radial 
displacement. The high-end tail of the BSL distribution shows a some- 
what opposite behaviour and its contribution to the surface cratering 
is eventually completely masked by the one of the ISL impactors. 



asteroids up to sizes of ~ 100 km (see figs. [TT1 and [T2"j) . BSL- related craters, whose 
sizes are generally distributed between 50 km and 150 km, are responsible for the 
high-end tail (150 — 250 km) of the simulated crater distributions (see figs. [TTland[T2l 
taking into account that the ISL-related crater frequencies are divided by 100). More 
in detail, for both asteroid we observe the following interesting result. In Jupiter's 
no- migration scenario the high-end tail (d > 100 km) of the crater size distribution 
is completely due to BSL-related impacts (see figs. [TT]and[T2l top left panels). For 
increasing values of the displacement of Jupiter's formation region from its present 
position, we expect only few BSL-related craters with size greater than 100 km (see 
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Figure 12. Normalised frequency versus crater size histograms of 
Ceres' impactors: red bars are those related to ISL impactors, blue 
ones to BSL impactors. The frequency of ISL impactors in all plots 
has been divided by a factor 100 to enhance the readability of the 
plots. Similarly to what we observed for Vesta in fig. [TTJ the high-end 
tail of the ISL crater size distribution moves towards bigger sizes for 
increasing values of Jupiter's radial displacement. The contribution of 
the BSL impactors to the surface cratering is eventually masked by the 
one of the ISL impactors, with the only exception of stochastic events 
(see bottom right panel). 



figs. [TTlCE^i top right and bottom left panels). At the same time, the high-end tail 
of the ISL-related crater distribution partially overlaps with the BSL-related crater 
distribution (ibid). In the 1AU migration scenario, the high-end tail of the ISL- 
related crater distribution completely overlap with the BSL-related one, cancelling 
its contribution in terms of crater frequency (see figs. [11] and [12j bottom right 
panels). In Ceres' case, we observe that in this migration scenario the contribution 
of BSL impacts could survive only for crater sizes higher than 200 km (see fig. [T21 
bottom right panel). 
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Figure 13. Normalised frequency versus impact energy (in units of 
autogravitation energy) histograms of Vesta's impactors: red bars are 
those related to ISL impactors, blue ones to BSL impactors. The 
frequency of ISL impactors in all plots has been divided by a factor 
100 to enhance the readability of the plots. The most energetic impacts 
(of the order of 10" 1 ) are due to BSL bodies in those scenarios where 
Jupiter's radial displacement is limited (no migration or 0.25 AU). ISL 
impacts in all scenarios and BSL impacts in the 0.5 AU and the 1 AU 
scenarios are limited to impact energies inferior to 10~ 2 . 

3.3.3. Impact energy. To evaluate the effects of ISL and BSL impacts on the internal 
structure of the asteroids, we computed the ratio between the energy delivered by 
the impacts and the self-gravitation energies of Vesta and Ceres (see figs. [13] and 
[TJJ. In Vesta's case the delivered energy spans in the range 10 -8 — 10 _1 (fig. [TBI 
while in Ceres' case, due to the greater size of the latter, the range is shifted by 
an order of magnitude, i.e. between 10~ 9 — 10~ 2 (fig. [Lit) . BSL impactors near 
the frequency peak in their distribution generally deliver 2 — 3 orders of magnitude 
more energy than ISL impactors at their corresponding peak (see figs. [TBI and IT%|) . 
In analogy to what we observed for the surface cratering, in Jupiter's no-migration 
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Figure 14. Normalised frequency versus impact energy (in units of 
autogravitation energy) histograms of Ceres' impactors: red bars are 
those related to ISL impactors, blue ones to BSL impactors. The 
frequency of ISL impactors in all plots has been divided by a factor 100 
to enhance the readability of the plots. Due to the bigger size of Ceres, 
the impact energies delivered by both ISL and BSL planetesimals in 
all migration scenarios are of the order of 10~ 2 or inferior. 



scenario ISL impacts dominate the low-energy tail of the distribution while BSL 
impacts are mainly responsible of the high-energy tail (ibid). For increasing values 
of Jupiter's radial displacement, however, the high-end tail of the ISL distribution 
shifts toward higher energy values, incresingly overlapping to the BSL distribution 
(ibid). The contribution of BSL impactors thus become more and more negligible 
for increasing values of Jupiter's radial displacement and is completely masked by 
that of ISL impactors in the 1AU migration scenario (ibid). This effect hold true 
also in Ceres' case, where the BSL population of impactors is enhanced due to the 
excitation of a mean motion resonance with Jupiter. From our results it is clearly 
visible as, due to its smaller size, impacts would deliver on Vesta an higher amount 
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of energy per units of self-gravitation energy (the high-energy tail being of the order 
of 10 _1 , see fig. [T3|) and are more likely to shatter or significally fracture the asteroid. 



4. Discussion 

Our results suggest that the cratering histories of Vesta and Ceres should be 
strongly correlated to the timescale of Jupiter formation and to its radial displace- 
ment. If Vesta or Ceres formed before Jupiter's core reached its critical mass, the 
expected crater distribution on the surfaces of the two asteroids should show differ- 
ent effects at different times. At the beginning, the crater distribution is dominated 
by impacts of small rocky bodies that should have lead to an intense and uniform 
craterisation. These low-velocity impacts could have contributed to the accretion of 
the asteroids. After Jupiter starts to accrete gas, a second phase of more energetic 
impacts takes place. Due to their lower number, BSL planetesimals contribute less 
to the craterisation but, with their higher impact energies, they affect the resurfac- 
ing and the thermal histories of the asteroids. A significant fraction of resonant ISL 
bodies should have impacted the asteroids at supersonic speed, contributing to their 
excavation and fragmentation. Jupiter's migration should have important effects on 
this picture, changing the impact speeds and the relative abundances of ISL and BSL 
populations. The marked dynamical differences between the families of impactors 
and the strong effect of Jupiter migration should help in interpreting and dating the 
different features that Dawn will observe on Vesta's and Ceres' surfaces. 
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